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I.  Implementation  of  a  Kinetic  Nucleation  Model  for  DSMC 

Comparison  of  First-time  Cluster  Nucleation 
Models  for  DSMC  Kinetic  Nucleation  vs  CNT 

Normalized  Cluster  number  density  Normalized  Cluster  Size 


Figure.  1 


Cluster  numbers  and  size  distributions  are  shown  for  a  simple  expansion  of  an 
Argon  gas  into  a  vacuum  through  a  sonic  nozzle  with  a  diameter  of  3.2  mm,  a  stagnation 
temperature  of  280  K,  and  a  stagnation  pressure  of  125  torr.  The  clusters  are  formed 
during  the  non-equilibrium  expansion  through  the  process  of  homogeneous  nucleation. 
Two  sets  of  comparisons  are  presented  for  cluster  number  density  and  average  cluster 
size,  with  the  top  contours  showing  the  kinetic  nucleation  result  and  the  bottom  contours 
showing  the  results  obtained  from  Classical  Nucleation  Theory  (CNT).  In  the  first 
comparison  of  cluster  number  density,  the  Kinetic  and  CNT  results  are  normalized  to  8 
xlO20  and  8  x  1017  m'3,  respectively.  In  the  comparison  of  the  cluster  size,  on  the  right 
hand  side,  the  contour  levels  are  normalized  to  50  and  1000  for  the  kinetic  and  CNT 
models,  respectively. 

The  most  important  conclusion  from  the  viewgraph  is  that  there  is  a  significant 
change  in  the  physics  of  a  homogeneous  condensation  flow  depending  on  which  cluster 
nucleation  model  is  used.  The  kinetic  nucleation  process  predicts  three  to  four  orders  of 
magnitude  more  clusters  are  formed  than  CNT.  Also,  the  kinetic  nucleation  model 
predicts  that  most  of  the  clusters  will  be  smaller  than  if  the  CNT  is  used.  There  are  also 
important  computational  consequences  when  switching  from  the  CNT  to  the  kinetic 
nucleation  model.  The  cluster  size  and  number  density  contour  distributions  are 
smoother  for  the  kinetic  nucleation  calculations  even  though  both  DSMC  calculations 
have  the  same  numerical  parameters  -  number  of  gas  and  cluster  particle  F„,  number  of 
cells  and  computational  domain.  The  shakiness  of  the  CNT  contours  is  due  to  the 
necessity  of  obtaining  the  local  gas  temperature  in  order  to  calculate  the  nucleation  rate. 
This  is  an  inherently  noisy  process  in  a  DSMC  calculation. 

Finally,  these  results  are  presented  for  Argon  because  there  is  a  data  base  of 
laboratory  measurements  to  which  we  can  validate  these  results.  This  activity  is  in 
progress.  When  the  method  is  extended  to  more  realistic  propellant  mixtures  of  water 
and  inorganic  impurities,  the  amount  of  nucleation  and  the  size  of  the  clusters  will  affect 
plume  signatures  and  the  degree  of  backflow  contamination. 

2.  Background  explanations: 

a.  DSMC  Implementation  of  cluster  formation:  A  two-dimensional, 
axisymmetric  DSMC  simulation  begins  from  the  imposed  starting  surface,  located  at 
approximately  one  orifice  diameter  downstream  of  the  orifice,  and  is  composed  of  45  jets. 
The  entire  DSMC  computational  domain  is  4  by  10  orifice  diameters  squared,  with  100  x 
250  cells,  each  cell  can  be  divided  up  to  16  subcells.  The  simulated  vapor  and  cluster 


particles  represent  1.0  x  10+,°  and  5.0  x  10+7  real  particles  respectively,  and  a  total  of 
about  860,000  simulated  particles  were  used.  A  time  step  of  2.0  x  10'9  sec  is  used  with 
about  24,000  steps  to  reach  the  steady  state  and  the  simulation  results  are  sampled  for 
200,000  steps. 

b.  Failure  of  CNT: 


The  DSMC/CNT  has  been  successful  in  predicting  average  cluster  properties  for  flows 
that  exhibit  very  low  amount  of  nucleation.  However,  it  fails  to  predict  cluster 
distribution  properties.  Since  homogeneous  condensation  had  never  been  incorporated 
directly  into  DSMC  it  was  still  a  major  advance  to  have  a  flow  coupled  condensation 
model.  However,  there  are  important  problems  with  the  inherent  CNT  assumptions  when 
modeling,  supersonic,  expanding  flows.  CNT  assumes  that  the  initial  clusters  are  bom  in 
a  constant  vapor  environment  and  can  only  grow  to  a  maximum  size.  The  initial  clusters, 
known  as  the  critical  clusters,  are  assumed  to  be  bom  in  an  equilibrium  isothermal  vapor 
environment.  The  classical  nucleation  rate,  J,  is  a  function  of  macroscopic  parameters 
as: 


V®*3  pw 


exp( 
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where  a,  the  cluster  surface  tension  and  pw,  the  cluster  mass  density,  are  the  most 
problematic  terms  for  small  clusters. 

c.  Implementation  of  a  full  kinetic  nucleation  model:  The  molecular  dynamics 
(MD)  approach  was  used  to  simulate  a  1-D  Argon  condensation  plume  along  the 
centerline.  A  geometry  criteria  and  history  tracking  method  were  used  to  analyze  how 
the  clusters  are  created  in  the  plume  and  it  was  validated  that  the  initial  dimers  were 
created  through  triple  collision  processes,  during  which  the  third  monomer  takes  away  the 
additional  energy.  The  triple  collision  can  be  expressed  as, 

A  +  A+A  — ^  >A2+A 


A  and  A2  are  atoms  and  stable  dimers  and 


where  P,  is  the  probability  of  creating  a  stable  dimer  from  a  triple  collision  and  /is  the 
energy  equilibrium  distribution  function. 


To  be  implemented  into  a  DSMC  scheme,  the  triple  collision  is  decomposed  into  two 
successive  binary  collision  steps: 

A  +  A  -»  A\ 


A2  +  A  — >  A2  +  A 

where  A 2  is  a  temporary  collision  complex  with  a  lifetime.  Complex-monomer 
collisions  are  referred  as  the  triple  collisions.  In  the  full  DSMC  simulation,  we  need  both 
P,  and  the  collision  complex  lifetime.  A  Hybrid  DSMC/MD  technique  was  employed  to 
obtain  P,  in  the  full  DSMC  simulations.  In  addition,  the  acceptance-rejection  principle  is 
used  to  determine  whether  a  binary  collision  occurs.  Each  successful  binary  collision  =>  a 
temporary  complex.  During  the  complex  lifetime,  complex-atom  pairs  are  selected  and 
the  acceptance-rejection  principle  is  applied  assuming  a  hard-sphere  like  interaction 
between  clusters  and  Argon  atoms.  The  collision  outcome  is  evaluated  according  to  the 
probability,  P,.  and  a  “successful”  complex-atom  collision  generates  a  stable  dimer. 
Subsequent  cluster  evolution  processes,  starting  from  dimers,  are  calculated  using  the 
previously  developed  models  for  cluster  sticking  and  evaporation. 

II:  Use  Of  Quasiclassical  Trajectory  Methods  To  Model  OH  Production 
Mechanisms  In  A  Side-Jet  RCS  Flow 


Use  of  Quasiclassical  Trajectory  Methods  to 
Model  OH  Production  Mechanisms  in  a  Side- Jet 
RCS  Flow 
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Figure  2. 


1.  Explanation  of  viewgraph: 


The  objective  of  this  work  is  to  understand  why  the  common,  readily  available 
Total  Collisional  Energy  (TCE)  model  for  chemical  reactions  in  DSMC  fails  for  the 
exchange  reaction  of  O+HCl  — »  OH+C1.  This  is  an  important  chemical  reaction 
occurring  between  the  efluents  of  a  solid-non-aluminized  reaction  control  systems  (RCS) 
and  free  stream  atomic  oxygen  present  in  large  quantities.  Initial  investigations  of  this 
multi-scale,  chemically  reacting  flow  showed  that  using  the  baseline  TCE  model 
approximately  10%  of  the  reaction  probabilities  would  be  unphysically  greater  than  unity, 
if  an  artificial  limiter  was  not  utilized  in  the  calculations.  For  this  reason  we  propose  to 
develop  a  new  chemical  model  that  will  be  valid  for  chemically  reactive,  hypervelocity 
collisions. 

We  use  the  Quasi-Classical  Trajectory/  Molecular  Dynamics  method  (QCT/MD), 
an  accurate  chemical  modeling  tool,  to  evaluate  the  reaction  and  total  collision  cross 
sections.  In  this  manner  we  can  determine  whether  it  is  possible  to  propose  a  correction 
to  the  TCE  model  and  we  investigate  the  sensitivity  of  the  full  DSMC  simulation  to  the 
different  chemical  models  for  this  reaction  for  a  hypersonic  condition  of  120  km  altitude, 
5  km/s,  a  typical  free  stream  condition  where  these  flows  occur. 


Background  explanations: 


The  classical  Hamiltonian  is  H=T+V=E,  where  T,  and  V  are  the  kinetic  and  potential 
energy  of  the  molecular  system  and  Hamilton’s  equations  of  motion  are 


dv  . 


P± 

m, 


where  i=l  ,2,3  (atoms),  j=l  ,2,3  (x,y,z);  py  and  ry  are  the  Cartesian  momenta  and 
coordinates,  and  m*  is  the  mass  of  the  i-th  atom.  We  solve  the  1 8  coupled  equations  for 
Pij  and  ry  using  a  4th-order  Runge-Kutta  method.  To  construct  the  initial  conditions  for 
each  trajectory  (i.e.,  internal  energy  of  target  HC1,  relative  velocity  of  O+HCl,  offset  or 
impact  parameter),  we  use  the  micro-canonical  sampling  to  obtain  the  initial  conditions  of 
the  scattering  event: 

W  =  Vfo,  -  V{R)\Pu  =  *>7  =  1.3 


where  V(R)  is  the  potential  energy  surface  (PES)  shown  on  the  viewgraph. 

To  obtain  final  state  distributions  of  product,  we  assume  that  the  vibrational  and 


rotational  energies  may  be  treated  classically  and  separable.  Each  trajectory  outcome  is 
analyzed  with  the  following  classical  relationships  for  angular  momentum  and  vibration. 


The  conclusions  of  this  work  is  that  both  the  reaction  and  total  cross  sections,  Cf, 
and  Or  ,  for  O+HCl  obtained  using  the  MD/QCT  method  affected  the  OH  production. 
Furthermore  it  was  found  that  Or  was  independent  of  the  HC1  internal  energy  and  the 
usual  Variable  Hard  Sphere  (VHS)  model  used  in  DSMC  could  be  corrected  for  high- 
temperature  O-HCl  collision  pairs  if  VHS  parameters  of  w=0.39,  drefi'i.9  A  at  1000  K  are 
used.  The  most  accurate  reaction  probability  calculated  by  the  QCT/MD  method, 
/V(MD),  was  less  than  0.4  at  120  km  altitude  for  freestream  velocity  of  5km/s. 

The  use  of  the  MD/QCT  model  is  limited  by  the  accuracy  of  the  PES  and  it  is  not 
typical  that  such  high  fidelity  PESs  are  available.  For  this  reason,  the  MD/QCT 
calculations  were  repeated  for  the  older,  more  readily  available,  but  less  accurate,  PES 
known  as  a  LEPS.  These  surfaces  could  be  more  easily  obtained  for  other  propel  lants- 
atomic  oxygen  collisions.  We  are  presently  analyzing  the  results  and  how  they  affect  the 
total  OH  production  in  the  complete  DSMC  simulation. 
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Abstract 

The  direct  simulation  Monte  Carlo  (DSMC) 
method  has  recently  been  developed  to  sim¬ 
ulate  homogeneous  condensation  in  a  free 
expansion  rocket  plume.  However,  cluster- 
monomer  and  cluster-cluster  collision  models 
as  well  as  the  determination  of  cluster  size 
were  simplified  in  the  previous  work  which  may 
greatly  impact  the  accuracy  of  numerical  sim¬ 
ulation  results.  In  this  work,  the  molecular 
dynamics  (MD)  method  is  used  to  simulate 
collision  and  sticking  probability  for  Ar  clus¬ 
ters.  The  Ashrigz-Poo  model  is  used  to  predict 
the  outcome  of  cluster-cluster  collisions.  These 
improved  models  are  then  integrated  into  a 
DSMC  code  to  predict  the  Rayleigh  scattering 
intensity  in  a  free  expanding  Argon  condensa¬ 
tion  plume  and  numerical  results  are  compared 
with  experimental  data  along  the  plume  cen¬ 
terline. 

1  Introduction 

Homogeneous  condensation  [1],  observed  in 
various  types  of  plumes  expanding  into  a  vac- 
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uum,  may  result  in  contamination  of  spacecraft 
surfaces  as  well  as  affect  the  observance  of  high 
altitude  plumes  [2]  [3].  Extensive  experimental 
studies  on  the  condensation  in  real  thrusters 
and  small  free-expanding  jets  have  been  con¬ 
ducted  since  the  1970’s  [4]  [5].  However,  it  is 
difficult  to  measure  cluster  size  and  number 
density  distributions  in  an  operational  plume. 
Therefore  a  prediction  tool  which  can  be  val¬ 
idated  with  small  thruster  data  will  provide 
crucial  information  about  homogeneous  con¬ 
densation  in  typical  chemical  rocket  plumes. 

Most  of  the  numerical  work  on  condensation 
coupled  to  a  gas  expansion  models  the  pro¬ 
cesses  in  ground-based  facilities  and  uses  con¬ 
tinuum  approaches  [6].  Since  free-expanding 
jets  are  mostly  in  the  transitional  to  rarefied 
regimes,  a  continuum  approach  would  be  in¬ 
applicable  to  simulate  condensation  in  a  free 
expansion  plume.  Thus,  a  more  appropriate 
kinetic  approach  should  be  used.  Such  an  ap¬ 
proach  must  address  the  processes  of  cluster 
nucleation,  growth,  decay  and  collisions  with 
other  clusters  and  monomers  as  well  as  the 
usual  gas  kinetics  and  is,  therefore,  computa¬ 
tionally  challenging. 

Recently,  the  direct  simulation  Monte  Carlo 
(DSMC)  [7]  method  has  been  applied  to 
simulate  homogeneous  condensation  in  free- 
expansion  plumes  [8]  [9].  The  previous  work  [8] 
can  be  summarized  as  follows.  Classical  nu¬ 
cleation  theory  (CNT)  [1]  was  used  to  model 
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the  generation  of  initial  nuclei,  which  are  as¬ 
sumed  to  be  created  at  the  local  critical  con¬ 
ditions.  Microscopic  models,  consistent  with 
CNT  theory,  have  been  developed  from  gen¬ 
eral  conservation  equations  to  model  sticking 
and  non-sticking  collisions,  as  well  as  evapora¬ 
tion  processes.  These  models  were  then  inte¬ 
grated  into  the  DSMC  method  to  simulate  con¬ 
densation  coupled  flow.  The  developed  DSMC 
condensation  model  was  numerically  validated 
by  comparing  the  simulation  results  with  ana¬ 
lytical  solutions  in  one  dimensional  test  cases, 
further  validated  by  comparison  with  Hagana’s 
scaling  laws  [4],  and  finally  applied  to  predict 
water  homogeneous  condensation  in  a  rocket 
plume. 

In  the  previous  work  [8],  clusters  were  mod¬ 
eled  as  hard  spheres,  with  radius  proportional 
to  is ,  where  i  is  the  number  of  molecules  in  the 
cluster.  The  probability  of  a  sticking  collision 
between  a  cluster  and  a  condensible  monomer 
was  simply  assumed  to  be  unity  [8]  or  0.1  [9] 
without  validation.  Cluster-cluster  collisions 
were  assumed  to  be  elastic  without  any  al¬ 
lowance  for  coalescence,  although  it  was  shown 
that  for  cases  considered  in  Ref.  [8]  the  number 
of  cluster-cluster  collisions  was  small.  Thus, 
there  is  room  for  improvement  of  the  previ¬ 
ous  numerical  simulations.  The  purpose  of  this 
work  is  to  improve  the  fidelity  of  the  models 
used  in  the  condensation  flow  by  use  of  the 
molecular  dynamics  (MD)  method. 

The  paper  is  organized  as  follows.  To  un¬ 
derstand  basic  cluster  model  and  interactions 
between  a  cluster  and  a  molecule,  we  choose 
a  molecular  dynamics  [10]  (MD)  approach  to 
simulate  cluster-molecule  collision  and  sticking 
processes  in  Sec.  2.  In  order  to  derive  gen¬ 
eral  conclusions  from  MD  simulations,  a  6-12 
Lennard- Jones  potential,  applicable  to  many 
species,  is  used  to  simulate  collision  processes 
between  Argon  clusters  and  Argon  monomers. 
Using  reduced  units,  the  MD  simulation  results 
of  Argon  species  are  compared  with  those  of 
nickel  from  Ref.  [11].  In  Sec.  3,  we  briefly 


review  the  microscopic  models  developed  in 
our  previous  work  [8]  and  explain  the  simi¬ 
larities  and  differences  of  the  previous  model 
with  the  corrected  models.  These  models  are 
then  used  in  Sec.  4  to  simulate  a  conden¬ 
sation  plume  in  the  direct  simulation  Monte 
Carlo  (DSMC)  method  of  a  free  expanding 
Argon  plume  through  a  sonic  orifice.  Us¬ 
ing  the  DSMC  results,  the  Rayleigh  scatter¬ 
ing  intensity  is  then  calculated  and  compared 
with  experimental  data  along  the  plume  cen¬ 
terline  [15]. 


2  Molecular  Dynamics 
Simulation 


2.1  Methodology 

Molecular  dynamics  (MD)  is  a  computer  sim¬ 
ulation  technique  which  allows  one  to  predict 
the  time  evolution  of  a  system  of  interacting 
particles,  e.g.,  atoms  or  molecules.  Each  real 
atom  in  MD  is  usually  modeled  as  a  point-size 
particle,  and  the  interaction  between  them  is 
modeled  by  an  interaction  potential.  Once  the 
system  initial  conditions  are  defined  in  terms 
of  particle  coordinates,  and  velocities,  the  sys¬ 
tem  evolution  in  time  is  followed  by  solving  the 
Hamilton  equations. 

The  traditional  6-12  Lennard-Jones  poten¬ 
tial,  valid  for  interactions  among  many  species, 
is  chosen  as  the  interaction  potential: 


V  =  4e 


er-ey 


(i) 


where  the  potential  parameters  for  Argon  are 
e  =  0.0103  eV  and  a  =  3.405  A  .  Note  that  in 
this  section,  we  specifically  focus  on  the  inter¬ 
actions  among  Argon  atoms. 

Trajectories  are  integrated  numerically  using 
the  Nordsieck  fifth-order  predictor-corrector 
algorithms  [10],  which  can  be  summarized  as 
follows.  First,  new  positions  as  well  as  time 
derivative  terms  up  to  fifth  order  are  predicted 
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at  a  new  time  step  using  a  Taylor  series  expan¬ 
sion.  Then,  the  force  acting  on  each  particle 
is  computed  for  the  predicted  positions,  and 
error  terms  are  calculated  by  comparing  force 
terms  between  the  previous  and  current  time 
steps.  Finally,  the  above  error  terms  are  used 
to  correct  the  particle  positions  and  derivatives 
for  the  next  time  step. 

2.2  Initial  Clusters 

To  simulate  cluster  collision  outcomes  using 
a  molecular  dynamics  approach,  the  initial 
conditions  (coordinates  and  velocities)  of  each 
molecule  within  a  cluster  must  be  specified  be¬ 
fore  the  simulation.  The  condition  of  dimers, 
obtained  from  previous  modeling  of  a  free  ex¬ 
pansion  plume  [12],  will  be  used  to  characterize 
the  initial  clusters  for  the  MD  simulations  dis¬ 
cussed  in  the  subsections  2.3  and  2.4. 

To  explore  the  physics  of  the  cluster  for¬ 
mation  processes  in  an  expanding  plume,  a 
molecular  dynamics  approach  was  used  to  sim¬ 
ulate  a  free  argon  expansion  along  the  core 
plume.  One  challenge  in  the  MD  modeling 
of  such  a  system  is  the  selection  of  a  crite¬ 
rion  to  distinguish  clusters,  especially  dimers, 
from  atoms  that  are  just  collision  pairs.  The 
Stillinger  criterion  [13]  and  a  history  tracking 
approach  were  used  to  find  stable  clusters  as 
follows.  First,  the  collision  duration  time  be¬ 
tween  two  atoms  were  simulated  by  the  MD 
approach  for  conditions  of  various  impact  pa¬ 
rameters  and  relative  velocities,  and  the  MD 
results  were  compared  with  Bunker’s  approxi¬ 
mate  formula[14].  It  was  found  that  the  molec¬ 
ular  collision  duration  time  is  usually  less  than 
lOps.  Then  each  molecular  position  in  the 
plume  was  output  every  lOps  during  the  MD 
simulation,  and  a  pure  geometric  Stillinger  cri¬ 
terion  was  used  to  identify  the  clusters  in  the 
plume.  In  this  way,  the  observed  clusters  were 
tracked  using  the  MD  output,  i.e.,  every  lOps, 
until  the  atoms  in  the  clusters  had  separated 
into  monomers.  We  found  that  the  dimer’s  life¬ 


time  could  last  more  than  30ps  in  the  plume, 
and  they  usually  originated  from  a  tri-atom 
collision  process  during  which  one  of  the  atoms 
takes  away  extra  energy,  while  the  other  two 
form  a  stable  dimer. 

To  evaluate  the  cluster  cross  section  and 
cluster-monomer  sticking  probabilities,  we  also 
need  to  prepare  clusters  of  size  greater  than 
two  for  cluster  configurations,  typical  of  ex¬ 
panding  flows.  Trimers  were  created  through 
a  dimer-monomer  sticking  collision,  tetramers 
through  trimer-monomer  sticking  collisions, 
and  so  on.  To  efficiently  generate  larger  clus¬ 
ters,  cluster-cluster  sticking  collisions  may  be 
used. 

Molecular  dynamics  simulation  results  of 
cluster-monomer  collisions  will  be  discussed  in 
the  following  subsection  2.3  and  2.4.  Each  MD 
simulation  case  is  repeated  for  various  trajec¬ 
tories  by  randomly  rotating  the  initial  cluster 
around  its  center  of  mass  and  the  results  av¬ 
eraged  over  the  trajectories  are  shown  in  the 
subsections  2.3  and  2.4. 

2.3  Cluster  Model 

The  cluster-monomer  collision  process  may  be 
used  to  define  the  cluster  cross  section,  as 
the  monomer-monomer  collision  defines  molec¬ 
ular  cross  section  in  the  traditional  DSMC  ap¬ 
proach.  To  characterize  the  cluster  model, 
the  collision  processes  of  cluster-monomers  are 
simulated  with  the  MD  method,  under  the  con¬ 
dition  of  various  impact  parameters,  b,  and 
relative  velocities,  vrei ■  Sticking  collisions  are 
obviously  counted  as  effective  collisions,  while 
for  non-sticking  collisions,  only  trajectories  in 
which  the  deflection  angle  is  at  least  10  de¬ 
grees  are  counted  as  collisions.  Note  that  the 
colliding  molecules  are  initially  separated  suffi¬ 
ciently  far  from  the  target  cluster  such  that  the 
effective  interaction  potential  between  them  is 
zero. 

The  collision  probability  for  each  case  (6, 
Vrei)  is  defined  as  the  ratio  of  the  number  of 


3 


trajectories  that  resulted  in  a  collision  to  the 
total  number  of  trajectories.  For  each  case,  a 
sample  of  200  trajectories  is  used  with  the  ori¬ 
entation  of  the  target  cluster  to  the  monomer 
randomly  chosen.  The  impact  parameter  b  in¬ 
creases  from  a  head-on  collision  value  of  zero  to 
a  value  in  which  the  cluster-monomer  collision 
probability  is  close  to  zero.  Relative  collision 
velocities  from  approximately  150  to  225  m/s 
are  chosen,  and  the  cluster  sizes  selected  for 
study  in  the  MD  simulations  ranges  from  a 
dimer  up  to  482-mers.  The  range  of  cluster 
size  and  collision  impact  parameter  are  typi¬ 
cal  in  a  condensation  flow  [8],  and  the  collision 
relative  velocity  corresponds  to  a  translational 
thermal  temperature  which  is  also  representa¬ 
tive  of  a  typical  nondimensional  temperature, 
the  ratio  of  the  collision  thermal  temperature 
to  the  species  freezing  temperature. 

Representative  MD  simulation  results  are 
shown  in  Fig.  1.  It  can  be  seen  that  the 
collision  probability  is  close  to  unity  for  the 
cases  with  relatively  small  impact  parameter. 
For  the  cluster  collision  probability  of  dimer- 
monomer  collisions,  for  example,  it  can  be  seen 
that  for  impact  parameters  b  less  than  bs  ~  7.5 
A  ,  the  probability  is  almost  independent  of 
the  impact  parameter  b.  For  collisions  with 
b  greater  than  bs,  the  collision  probability  de¬ 
creases  rapidly  to  zero  as  the  impact  parameter 
increases.  It  can  also  be  seen  that  bs,  defined 
as  the  critical  impact  parameter,  increases  as 
the  cluster  size  increases. 

The  maximum  impact  parameter  bm,  cho¬ 
sen  here  to  correspond  to  a  collision  probabil¬ 
ity  of  0.1,  can  be  used  to  calculate  the  cluster- 
monomer  collision  cross  section,  a,  as  a  =  tv b^. 
It  is  related  to  the  cluster  and  monomer  diam¬ 
eters  as  [7] 

u  _dc  +  d0  /<v, 

bm  —  2  (2) 

where  dc  and  d0  are  the  cluster  and  molecule  di¬ 
ameters,  respectively.  Thus,  the  cluster  radius, 
Ri  =  can  be  determined  from  the  MD  sim¬ 
ulation  results,  since  the  monomer  diameter  is 


known.  The  MD  results  are  shown  in  Fig.  2 
with  a  relative  collision  velocity  of  182 m/s  and 
are  fitted  to  an  analytic  general  cluster  radius 
equation  [16]: 

Ri=Ai*+B  (3) 

where  i  is  the  number  of  molecules  in  the  clus¬ 
ter,  and  A  and  B  are  constants.  It  can  be  seen 
that  the  MD  simulation  results  for  Argon  clus¬ 
ters  fit  this  equation  quite  well  by  choosing  A 
2.3A  and  B  3.4A  (dashed  line  in  Fig.  2). 

The  molecular  dynamics  simulations  also 
show  that  cluster-monomer  collision  probabil¬ 
ities  are  a  function  of  the  relative  collision 
velocities.  Results  of  collision  probabilities 
for  dilferent  cluster  sizes  and  relative  cluster 
monomer  velocities  are  shown  in  Fig.  3.  It 
can  be  seen  that  the  critical  impact  parame¬ 
ter  increases  as  the  relative  collision  velocity 
decreases,  and  therefore,  the  cluster  cross  sec¬ 
tion  would  also  increase  as  the  relative  velocity 
decreases.  To  model  the  dependence  of  clus¬ 
ter  diameter  on  the  relative  collision  velocity 
cv,  we  use  the  variable  hard  sphere  (VHS)  [7] 
model 

d  =  dTef(^±y  (4) 

Cr 

where  the  subscript  ref  refers  to  a  reference 
value  and  v  is  a  constant  related  to  species 
viscosity.  The  MD  results  for  cluster  radii,  cal¬ 
culated  for  various  relative  collision  velocities, 
are  shown  in  Fig.  4,  and  the  simulation  results 
fitted  to  Eq.  (4)  are  shown  as  dashed  lines.  Fig¬ 
ure  4  shows  that  the  cluster  radii  are  approx¬ 
imately  inversely  proportional  to  the  relative 
collision  velocity,  and  the  slopes  for  different 
cluster  sizes  are  similar.  The  reference  cluster 
diameter,  dref,  may  be  estimated  from  Eq.  (3). 
For  larger  cluster  size,  the  dref  is  higher  so  that 
a  smaller  v  is  required  to  correct  its  diameter 
for  the  impact  of  relative  collision  velocity.  It 
can  be  seen  that  as  the  cluster  size  increases 
from  dimer  to  147-mers,  the  constant  v  de¬ 
creases  approximately  from  0.3  to  0.07.  Note 
that  for  Argon  vapor,  the  constant  v  is  about 
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0.31.  Thus,  the  impact  of  collision  relative  ve¬ 
locity  on  the  cluster  cross  section  becomes  less 
as  the  cluster  size  increases,  which  indicates 
that  the  hard  sphere  model  may  be  used  to 
model  relatively  large  clusters  in  the  DSMC 
simulation.  As  Fig.  4  shows,  the  error  in  us¬ 
ing  the  hard  sphere  model  increases  for  small 
clusters. 

2.4  Sticking  Collision  Probabil¬ 
ity 

The  sticking  probability  for  a  collision  process 
between  a  cluster  and  a  molecule  can  be  simu¬ 
lated  in  a  manner  similar  to  the  collision  prob¬ 
ability.  Initially,  the  target  Argon  cluster  is  at 
the  origin  in  a  static  state1,  while  the  collid¬ 
ing  Argon  atom  is  approximately  50  A  away 
from  the  cluster.  It  takes  about  25  ps  before 
the  atom  collides  with  the  target  cluster.  The 
simulation  time  for  a  colliding  case  is  approx¬ 
imately  100  ps,  which  is  also  long  enough  to 
observe  whether  the  outcome  is  a  sticking  or 
non-sticking  collision.  To  get  good  statistics, 
we  perform  2500  trajectories  for  each  colliding 
case  by  randomly  rotating  the  target  cluster, 
as  discussed  in  Sec.  2.3. 

It  can  be  expected  that  cluster-monomer 
sticking  probability  is  also  dependent  on  the 
relative  collision  velocity,  especially  for  small 
clusters.  Similarly  to  the  collision  probability, 
the  impact  of  the  relative  velocity  on  the  stick¬ 
ing  probability  decreases  as  the  cluster  size  in¬ 
crease.  Figure  5  shows  the  MD  simulation  re¬ 
sults  of  sticking  probabilities  for  various  cluster 
sizes  and  relative  velocities  as  a  function  of  im¬ 
pact  parameter.  This  figure  indicates  that  the 
sensitivity  of  the  sticking  probability  decreases 
as  the  cluster  size  increases,  and  it  may  be 
neglected  for  the  clusters  larger  than  10-mers 
with  good  approximation.  To  highlight  the  size 
dependence  of  the  sticking  probability,  Fig.  6 

lA  static  state  indicates  one  with  no  translational 
velocity,  and  only  thermal  motion  or  internal  energy. 


shows  the  sticking  probability  as  a  function  of 
impact  parameter  for  different  size  clusters  all 
with  a  relative  collision  velocity  of  182 m/s.  It 
can  be  seen  that  the  sticking  collision  probabil¬ 
ity,  similar  to  the  collision  probability,  exhibits 
a  critical  impact  parameter  dependence.  For 
collisions  with  impact  parameters  less  than  the 
critical  impact  parameter,  the  sticking  proba¬ 
bilities  are  almost  constant,  and  are  unity  for 
cluster  sizes  larger  than  10-mers.  While  for  col¬ 
lisions  with  impact  parameters  larger  than  the 
critical  impact  parameter,  the  sticking  proba¬ 
bilities  decrease  to  zero  more  quickly  than  the 
collision  probabilities.  As  can  be  seen  from 
Fig.  6,  the  critical  collision  impact  parameter 
for  the  sticking  collision  probability  increases 
as  the  cluster  size  increases. 

In  the  DSMC  simulation  of  an  expanding 
flow,  consisting  of  gas  and  clusters,  each  col¬ 
lision  is  not  modeled  in  detail  and  only  aver¬ 
age  collision  properties  are  used.  For  this  rea¬ 
son,  the  dependence  of  the  sticking  probability 
on  impact  parameter  must  be  averaged  over 
the  collision  cross  section.  The  average  stick¬ 
ing  probability,  P*„,  can  be  calculated  from  the 
MD  results  for  each  cluster  based  on  a  b 2  dis¬ 
tribution  of  collision  pairs, 

fr»,m 

pL  =  db 2  (5) 

where  is  the  cluster  maximum  impact  pa¬ 
rameter  or  collision  radius,  given  by  Eq.  (2). 
Figure  7  shows  the  averaged  sticking  proba¬ 
bility,  Plav,  as  a  function  of  cluster  size  for  a 
relative  collision  velocity  of  182 m/s.  It  can  be 
seen  that  for  cluster  size  smaller  than  10,  the 
sticking  collision  probability  increases  quickly 
as  the  cluster  size  increases;  while  for  cluster 
size  larger  than  10,  the  sticking  collision  proba¬ 
bility  is  independent  of  cluster  size  and  is  very 
close  to  unity. 

Therefore,  the  MD  simulation  results  show 
that  for  clusters  larger  than  10-mers,  cluster- 
monomer  interactions  may  be  modeled  approx- 
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imately  with  the  hard  sphere  model,  with  a 
constant  sticking  coefficient  of  unity.  The  MD 
simulation  study  confirms  the  selection  of  these 
models  utilized  in  earlier  work  [8]. 


2.5  Comparison  of  MD  results 
for  different  Lennard- Jones 
Systems 

Collision  processes  between  argon  clusters  and 
molecules  have  been  simulated  in  the  previous 
sections,  and  important  results  of  cluster  size 
and  sticking  coefficients  are  shown  in  Fig.  2  and 
Fig.  7  respectively.  Since  Lennard-Jones  po¬ 
tential  is  generally  used  to  describe  molecular 
interaction  for  many  other  species,  the  simu¬ 
lated  results  for  Argon  may  be  able  to  roughly 
estimate  other  species  without  conducting  nu¬ 
merically  identical  simulations  using  various 
parameters.  To  validate  our  MD  simulations, 
we  would  like  to  compare  our  numerical  results 
with  other  independent  data  [11]. 

Nickel  cluster  collision  cross  sections  and 
sticking  collision  process  has  been  simulated 
by  the  MD  method  in  the  Ref.  [11],  using  a 
Lennard-Jones  potential.  First  we  may  di¬ 
rectly  compare  Ni  and  Ar  cluster  sticking  prob¬ 
abilities,  which  is  a  nondimensional  parameter. 
Note  that  the  nickel  cluster  sticking  probabil¬ 
ity  is  defined  as  the  ratio  of  the  cluster  stick¬ 
ing  cross  section  to  the  total  cross  section  and 
both  quantities  are  given  in  Ref.  [17].  The 
comparison  of  Argon  and  Nickel  cluster  stick¬ 
ing  probability  is  shown  in  Fig.  7,  represented 
by  open  circles  and  solid  triangles  respectively. 
Figure  7  suggests  that  the  trends  for  small  clus¬ 
ter  sticking  probability  discussed  in  Sec.  2.4  is 
similar  for  both  argon  and  nickel  clusters. 

To  compare  cluster  collision  cross  sections 
for  different  species,  we  have  to  use  reduced 
units  to  analyze  the  MD  simulation  results  of 
Ref.  [17].  We  assume  that  the  cluster  cross 


section  radius  can  be  expressed  as: 

ep 

rt  +  ri  oc  af(e,m)  <x  a — -  (6) 

where  r*  is  the  radius  of  a  cluster  consisting  of 
i  monomers,  r\  is  the  monomer’s  radius,  a  and 
€  are  the  Lennard-Jones  parameters,  m  is  the 
molecular  mass,  and  p  and  q  are  constants. 

The  nickel  cluster  cross  sections  are  simu¬ 
lated  in  Ref.  [11]  with  two  Lennard-Jones  po¬ 
tentials,  the  bulk  fitted  (LJB)  and  dimer  fit¬ 
ted  (LJD)  potential.  The  simulation  results 
for  clusters  of  size  from  dimer  to  13-mers  are 
given  in  Ref.  [11],  e.g.  a  12-mers  cluster  cross 
section  is  about  250  A  2  for  the  LJB  potential, 
and  200  A  2  for  the  LJD  potential.  The  con¬ 
stant  p  in  Eq.  (6)  can  be  estimated  by  compar¬ 
ing  those  two  simulated  12-mers  cross  sections 
and  was  found  to  be  approximately  0.18.  Ex¬ 
amination  of  Fig.  2  shows  that,  the  12-mers  ar¬ 
gon  cluster  cross  section  is  approximately  363 
A  2  which  compared  with  the  Nickel  cluster 
computed  with  the  LJB  potential,  allows  one 
to  estimate  a  q ,  value  of  approximately  1.56. 

By  using  Eq.  (6),  the  Argon  cluster  colli¬ 
sion  cross  sections  can  be  used  to  predict  nickel 
cluster  collision  cross  sections.  A  comparison 
of  the  Ar  derived  predicted  Ni  cross  sections 
with  those  obtained  in  Ref.  [11]  is  shown  in 
Fig.  8.  It  can  be  seen  that  the  argon  cross  sec¬ 
tions  can  predict  nickel  cluster  cross  sections 
well  for  both  Ni  interaction  potentials.  Thus, 
the  Argon  cluster  sticking  probability  and  col¬ 
lision  cross  section  are  consistent  with  those 
of  nickel  clusters,  and  these  results  may  be  ex¬ 
tended  in  the  future  to  approximately  estimate 
other  species  described  by  theLennard-Jones 
potential,  such  as  water. 

3  Microscopic  Models 

Based  on  the  general  conservation  equations, 
the  developed  microscopic  models  had  been  in¬ 
tegrated  into  a  DSMC  computational  scheme 
to  predict  homogeneous  condensation  in  a 
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rocket  plume.  Detail  descriptions  of  these 
models  can  be  found  in  the  previous  work  [8]. 
To  predict  Argon  condensation  in  a  free  ex¬ 
panding  plume  in  Sec.  4,  we  briefly  review 
these  models,  and  discuss  the  improvements  in 
detail  as  follows. 


3.1  Simulated  Molecular  Models 

To  simulate  condensation  coupled  flow,  all 
species  including  monomers  and  clusters  are 
represented  by  simulated  molecules  in  the 
DSMC  approach.  The  variable  hard  sphere 
(VHS)  model  is  used  to  model  simulated 
molecules  representing  real  atoms  or  molecules 
in  the  plume,  and  the  Larsen-Borgnakke  model 
is  used  to  redistribute  molecule-molecule  colli¬ 
sion  energy,  based  on  the  number  of  degrees  of 
freedom.  Note  that  both  rotation  and  vibra¬ 
tion  energy  transfer  are  considered  in  molecule- 
molecule  collisions. 

An  analytic  cluster  hard  sphere  model  is  de¬ 
scribed  in  Ref.  [16],  providing  a  general  asymp¬ 
totical  equation  for  cluster  radius  similar  to 
Eq.  (3).  Briehl  et  al  [19]  applied  this  cluster 
model  to  simulate  small  Cu  cluster  growth  and 
decay  processes  using  the  DSMC  approach. 
This  cluster  model  has  been  validated  by  our 
MD  simulations  for  Argon  cluster,  discussed  in 
Sec.  2.2.  It  was  further  shown  that  Argon  clus¬ 
ter  model  is  consistent  with  other  independent 
studies  of  Nickel  clusters  [11]  in  Sec.  2.5.  Thus, 
we  will  use  the  hard  sphere  cluster  model  in  our 
DSMC  simulation. 


3.2  Nucleation  Model 

Following  CNT  theory,  initial  clusters,  called 
nuclei,  are  introduced  into  the  computational 
domain.  The  nuclei  have  the  same  size  as  the 
local  critical  clusters,  and  they  are  in  an  equi¬ 
librium  state  with  the  surrounding  gas.  The 
nucleation  rate,  J,  is  given  in  the  CNT  theory 


as: 


where  a  is  the  cluster  surface  tension,  m  is 
molecular  mass,  pv  is  vapor  density,  pw  is  clus¬ 
ter  density,  and  r*  is  the  local  critical  cluster 
radius.  Thus,  during  a  time  step  At  in  a  cell 
with  volume  AV,  the  number  of  new  simulated 
nuclei  Nc  is  calculated  as 

Nc  =  JcAVAt/Fc  (8) 

where  Fc  is  the  number  of  real  molecules 
represented  by  a  simulated  cluster  molecule. 
The  latent  heat,  generated  during  the  nucle¬ 
ation  process,  is  distributed  to  the  ambient  gas 
molecules  in  a  cell  according  to  energy  equipar- 
tition  principle.  Using  the  bulk  theory,  the  la¬ 
tent  heat,  El,  is  given  by 

El  =  LMi  (9) 

where  L  is  the  specific  latent  heat  of  evapora¬ 
tion,  and  Mi  is  the  mass  of  a  cluster  consisting 
of  i  monomers. 

3.3  Collision  Models 

The  collision  process  between  a  cluster  and  a 
foreign  molecule  is  referred  to  as  non-sticking 
collisions,  while  the  outcome  of  a  collision  be¬ 
tween  a  cluster  and  its  own  molecule  may 
be  sticking,  or  non-sticking,  depending  on  the 
sticking  collision  probability.  Based  on  molec¬ 
ular  dynamics  simulations,  a  general  expres¬ 
sion  for  both  argon  and  water  homogeneous 
sticking  coefficient,  q,  has  been  summarized 
as[21] 


where  V1  and  V9  are  the  specific  volumes  of  gas 
and  liquid  and  E0  is  the  energy  difference  of  the 
minimum  quantum  level  between  the  activated 
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complex  and  gas.  Equation  (10)  shows  that 
sticking  coefficient  is  close  to  unity  in  the  low 
temperature  region  and  decreases  as  the  tem¬ 
perature  increases.  As  shown  in  our  MD  simu¬ 
lations,  sticking  coefficient  for  clusters  smaller 
than  11-mers  may  not  agree  with  the  predic¬ 
tion  of  Eq.  (10)  since  it  is  based  on  the  con¬ 
dition  of  a  flat  liquid  surface.  The  impact  of 
correctly  modeling  the  sticking  probability  of 
clusters  smaller  than  11-mers  on  the  conden¬ 
sation  results  will  be  discussed  in  Sec.  4. 

In  the  previous  work,  cluster-cluster  colli¬ 
sions  are  assumed  to  be  elastic  collision,  with¬ 
out  considering  the  coalescence  effect.  In  this 
work,  the  semiempirical  Ashgriz-Poo  model  [22] 
is  chosen  to  predict  the  outcome  of  cluster- 
cluster  collisions  that  may  lead  to  coalescence 
or  separation.  The  Weber  number,  defined  as 
the  ratio  of  collision  kinetic  energy  to  droplet 
surface  energy,  provides  a  cluster-cluster  col¬ 
lision  separation-coalescence  boundary.  The 
Weber  number,  We,  is  given  by: 

We=P^A  (H) 

<?d 

where  Pd  is  the  cluster  density,  di  is  the  diam¬ 
eter  of  the  smaller  cluster,  Cr  is  the  cluster- 
cluster  relative  velocity,  and  04  is  the  clus¬ 
ter  surface  tension.  The  Ashgriz-Poo  model 
has  been  verified  for  water  droplets  by  MD 
simulations  using  a  6-12  Lennard- Jones  po¬ 
tential,  and  was  used  to  predict  the  outcome 
of  droplet  collisions  in  an  inductively  coupled 
plasma  [23].  The  detail  description  of  the 
Ashgriz-Poo  model  can  be  found  in  Ref.  [22]. 

4  Results  and  Discussions 

The  MD  simulation  results,  cluster  collision 
cross  sections  and  cluster-monomer  sticking  co¬ 
efficients,  will  be  integrated  into  the  DSMC 
simulations.  To  understand  the  impact  of 
these  improvements  and  check  the  simulation 
results  presented  in  the  previous  work  [8],  we 


want  to  first  compare  the  improved  DSMC 
simulation  results  with  previous  results. 

A  free  argon  expanding  flow  through  a  sonic 
nozzle  with  a  convergent  angle  of  30  degrees 
into  a  vacuum  environment  was  simulated  in 
Ref.  [8].  The  plume  has  a  stagnation  pres¬ 
sure  of  5.30  x  104  Pa,  and  a  stagnation  tem¬ 
perature  of  180  K.  The  DSMC  simulation  pa¬ 
rameters  and  techniques  for  this  case  were  de¬ 
scribed  in  detail  in  Ref.  [8].  Comparisons 
between  improved  (new)  and  previous  results 
(old)  are  shown  in  Fig.  9.  Since  the  MD  cor¬ 
rected  cluster-monomer  collision  cross  sections 
are  larger  than  the  previous  ones,  the  compar¬ 
ison  of  average  cluster  size,  shown  in  the  top 
of  the  figure,  suggests  that  the  cluster  size  ob¬ 
tained  by  the  improved  model  is  larger  than 
the  previous  size.  However,  after  the  cluster 
size  is  normalized  to  the  corresponding  termi¬ 
nal  cluster  size,  the  middle  figure  of  Fig.  9 
shows  that  the  normalized  curves  of  the  im¬ 
proved  and  previous  models  are  almost  the 
same.  The  middle  figure  suggests  that  the  im¬ 
proved  model  still  validates  the  scaling  laws 
of  Hagena  [4],  as  was  shown  in  the  previous 
work  [8].  Finally,  the  cluster  number  densities 
are  compared  in  the  bottom  of  Fig.  9  and  the 
results  suggest  that  the  cluster  number  density 
obtained  by  the  improved  model  is  slightly  less 
than  the  previous  value.  The  main  reason  of 
the  decreased  cluster  number  density  is  that 
small  clusters  in  the  improved  model  have  had 
more  chance  to  evaporate  into  separate  atoms 
due  to  relative  smaller  sticking  coefficients. 

Next,  the  improved  condensation  model 
is  applied  to  a  new  expanding  condensation 
plume  experimentally  studied  in  the  Arnold 
Engineering  Development  Center  (AEDC)  in 
the  1970’s  [15].  Extensive  sets  of  Rayleigh 
scattering  intensity  data  were  measured,  which 
may  be  directly  compared  with  numerical  re¬ 
sults.  A  brief  summary  of  the  scattering  inten¬ 
sity  is  given  below  and  detail  information  can 
be  found  in  Ref.  [15]. 

Assuming  that  the  condensing  flow  field 


8 


is  composed  of  a  collection  of  gas  phase 
monomers  and  i-mers  with  the  respective  num¬ 
ber  density  N{  and  polarizability  or*,  the  ratio 
of  scattered  intensity,  Is,  to  the  incident  inten¬ 
sity,  is  given  by: 

Is/I0  =  ^KNtayX40  (12) 

where  it'  is  a  coefficient  containing  transmis¬ 
sion  and  calibration  factors,  and  Ao  is  the  wave¬ 
length.  Note  that  polarizability  a*  for  i-mers 
is  assumed  to  be  i  times  larger  than  molecular 
polarizability  cq  [15]. 

Here  we  will  use  the  DSMC  method  to  simu¬ 
late  condensation  in  a  free  expanding  pure  ar¬ 
gon  plume  through  a  sonic  orifice,  which  had 
experimentally  done  in  AEDC  [15]  with  a  stag¬ 
nation  pressure  of  250  torr,  stagnation  temper¬ 
ature  of  280  K  and  orifice  diameter  of  3.2  mm. 
The  details  of  DSMC  simulation  models  and 
techniques  are  developed  in  Ref.  [8],  which  can 
be  briefly  summarized  as  follows.  The  origi¬ 
nal  2D  axisymmetric  SMILE  code  [24]  is  mod¬ 
ified  to  simulate  nucleation,  cluster-molecule 
and  cluster-cluster  collision,  and  evaporation 
processes.  The  simulation  domain  extends  4 
and  10  times  the  orifice  diameter  in  the  radial 
and  axial  directions,  with  100  and  250  cells  re¬ 
spectively.  Each  cell  can  be  divided  into  up  to 
16  subcells  according  to  the  flow  gradients  ob¬ 
tained  during  the  simulation,  and  a  time  step 
of  2.0  x  10~9  s  is  used.  To  decrease  the  compu¬ 
tational  cost,  radial  weights  are  used  for  dis¬ 
tributing  the  number  of  simulated  molecules 
evenly  in  the  radial  direction  for  monomers 
and  clusters.  A  cluster  weighting  factor  of 
1.0  x  10~5  is  used  in  the  simulation  due  to  the 
several  orders  of  magnitude  difference  in  the 
cluster  and  vapor  number  density.  Since  the 
non-condensation  region  close  to  the  orifice  is 
too  dense  to  be  modeled  by  DSMC,  the  simu¬ 
lation  begins  from  a  starting  surface  which  is 
created  with  the  continuum  solver  GASP  [25]. 
About  3,200,000  simulated  molecules  are  used 
to  represent  clusters  and  vapor  molecules  at 
steady  state.  A  typical  simulation  takes  72 


hours  on  18  parallel  AMD  Athlon  processors 
of  CPU  1526  MHz. 

The  DSMC  simulation  results  for  the  con¬ 
densation  plume  are  shown  in  Fig.  10,  includ¬ 
ing  vapor  number  density  contours  (top),  clus¬ 
ter  number  density  contours  (middle),  and  av¬ 
erage  cluster  size  contours  (bottom).  It  can  be 
seen  that  the  initial  cluster  or  nuclei  appear  in 
a  nucleation  region  where  cluster  number  den¬ 
sity  increases  quickly  along  the  flow  direction 
consistent  with  the  increased  degree  of  super¬ 
saturation  in  the  free  expanding  process.  Up¬ 
stream  of  this  region,  the  gas  is  essentially  an 
expanding  plume  without  condensation.  Fur¬ 
ther  downstream  of  the  nucleation  region,  clus¬ 
ter  number  density  decreases  mainly  due  to  the 
expansion.  Since  the  cluster  number  density 
is  only  on  the  order  of  1017  per  m3,  there  are 
few  cluster-cluster  collisions  so  that  the  impor¬ 
tance  of  the  coalescence  process  on  the  cluster 
number  density  is  limited.  The  average  cluster 
size  is  usually  below  10-mers  in  the  nucleation 
region,  while  it  increases  downstream  of  the 
nucleation  region  due  to  the  sticking  collision 
process  between  cluster  and  vapor  molecules 
and  cluster  evaporation  processes. 

The  distributions  of  average  cluster  size, 
cluster  and  vapor  number  density  along  the 
core  flow  are  shown  in  Fig.  11.  The  figure 
shows  that  vapor  number  density  is  about  5 
orders  of  magnitude  greater  than  the  clus¬ 
ter  number  density,  while  the  maximum  clus¬ 
ter  size  in  the  computational  domain  is  about 
300-mers.  The  average  cluster  size  seems  to 
grow  continuously  beyond  the  computational 
domain  until  the  vapor  number  density  be¬ 
comes  so  low  that  the  contribution  of  sticking 
collisions  to  cluster  growth  process  counteracts 
the  cluster  evaporation  effect.  The  Rayleigh 
scattering  intensity  is  calculated  along  the  core 
flow  using  Eq.  (12),  and  is  compared  with  ex¬ 
perimental  data  in  Fig.  12.  Note  that  the 
dashed  line  in  Fig.  12  indicates  the  change 
in  Rayleigh  scattering  intensity  along  the  core 
of  a  plume  for  which  condensation  is  ignored. 
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Upstream  of  the  nucleation  region,  there  is 
no  condensation  in  the  plume  which  leads  to 
a  Rayleigh  scattering  intensity  decrease  along 
the  flow  direction  due  to  rarefaction.  How¬ 
ever,  the  emergence  of  clusters  and  the  cluster 
growth  process  counteracts  the  decrease  of  va¬ 
por  number  density  in  the  condensation  region, 
leading  to  a  scattering  intensity  downstream 
of  the  nucleation  region  that  is  larger  than  its 
corresponding  value  in  the  non-condensation 
plume.  Figure  12  also  shows  that  our  numer¬ 
ical  results  agree  reasonably  well  with  the  ex¬ 
perimental  data,  thus,  we  expect  the  DSMC 
condensation  model  may  predict  the  conden¬ 
sation  phenomenon  in  a  free  expanding  plume 
correctly.  Note  that  the  appearance  of  nuclei 
in  the  DSMC  simulation  is  delayed  1.2  x  10-5 
sec  to  make  the  initiation  of  condensation  in 
the  simulation  coincide  with  the  onset  of  con¬ 
densation  at  the  position  indicated  in  the  ex¬ 
periment. 


Figure  13  compares  the  Rayleigh  scatter¬ 
ing  intensity  contours  between  a  condensation 
(top)  and  non-condensation  (bottom)  plume. 
For  the  non-condensation  case,  the  plume  scat¬ 
tering  intensity  decreases  monotonically  in  the 
expanding  directions,  due  to  the  decrease  in 
the  vapor  number  density.  However,  for  the 
condensation  plume,  the  scattering  intensity 
decreases  more  slowly  in  the  axial  direction 
after  the  initiation  of  condensation.  Consis¬ 
tent  with  the  results  shown  in  Fig.  11,  there  is 
a  counteracting  effect  between  the  increase  of 
the  cluster  size  and  the  decrease  of  the  cluster 
number  density.  Our  numerical  simulation  re¬ 
sults  show  that  there  would  be  a  long  tail  in 
the  Rayleigh  scattering  intensity  distribution 
in  the  free  expansion  condensation  flow.  The 
exact  spatial  distribution  should  be  further  val¬ 
idated  by  future  experiments. 


5  Summary 

The  general  molecular  dynamics  (MD)  ap¬ 
proach  is  used  to  verify  an  analytical  clus¬ 
ter  model,  and  cluster  sticking  collision  model. 
The  semiempirical  Ashriz-Poo  model  is  used  to 
predict  the  outcome  of  cluster-cluster  collision. 
The  use  of  LJ  potential  reduced  units  shows 
that  MD  simulation  results  of  Argon  cluster 
agree  well  with  those  of  Nickel  cluster  done  in 
an  other  independent  work,  indicating  that  the 
MD  results  in  this  work  may  be  approximately 
transferred  to  other  species. 

The  MD  simulation  results  are  applied  into 
the  DSMC  model,  which  is  developed  to  sim¬ 
ulate  homogeneous  condensation  in  a  free  ex¬ 
panding  Argon  plume  through  a  sonic  orifice. 
According  to  the  DSMC  simulation  results, 
Rayleigh  scattering  intensity  is  estimated  and 
compared  with  the  experimental  data  suggest¬ 
ing  that  the  condensation  model  developed  in 
this  work  may  reasonably  well  predict  homoge¬ 
neous  condensation  in  a  free  expanding  plume. 
The  Rayleigh  scattering  contours  in  a  conden¬ 
sation  plume  are  predicted,  which  may  need  to 
be  validated  by  experimental  data  in  the  fu¬ 
ture. 
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Figure  1:  Comparison  of  collision  probability  for  some  clusters  from  dimer  to  482-mers  under  the 
condition  of  various  impact  parameters,  with  a  relative  collision  velocity  of  182 m/s. 
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Figure  2:  Average  cluster  radius  from  MD  simulations  with  a  relative  collision  velocity  of  182 m/s 
is  fitted  to  an  approximated  Equation,  j Rj  =  2.3 i*  +  3.4. 
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Figure  3:  Comparison  of  collision  probability  for  some  clusters  under  the  condition  of  various 
impact  parameters  and  relative  collision  velocities. 
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Figure  4:  Cluster  radius  simulated  from  MD  are  fitted  with  equations:  r  =  rre/( 
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Figure  5:  MD  simulation  results  of  sticking  collision  coefficient  or  probabilities  for  some  clusters 
under  the  condition  of  various  impact  parameters  and  relative  collision  velocities. 
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Figure  6:  MD  simulation  results  of  sticking  collision  coefficient  for  clusters  of  different  sizes  as  a 
function  of  impact  parameter,  with  a  relative  collision  velocity  of  182m/s. 


Figure  7:  Average  cluster  sticking  probabilities,  for  argon  clusters  from  dimer  to  76-mers  with  a 
relative  collision  velocity  at  182m/s,  are  indicated  by  open  circles  compared  with  nickel  cluster 
sticking  probabilities  as  solid  symbols. 
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Cluster  Size 

Figure  8:  Comparison  of  predicted  nickel  cluster  cross  section  from  simulated  argon  clusters  using 
Eq.  (6)  with  MD  simulation  results. 
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Figure  9:  Contours  of  Ar  cluster  size  (top),  normalized  cluster  size  (middle),  and  cluster  number 
density  (bottom)  along  the  plume  centerline  between  the  improved  (new)  and  previous  (old) 
results.  Number  density  is  given  in  per  m3,  Axial  position  is  normalized  to  the  nozzle  throat 
radius  R. 
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Figure  10:  Contours  of  Ar  vapor  number  density  (top),  cluster  number  density  (middle)  and 
average  cluster  size  (bottom)  in  the  condensation  plume,  discussed  in  Sec.  4.  Number  density  is 


given  m  per  m  . 
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Figure  11:  The  distribution  of  Argon  gas  and  cluster  number  density,  and  average  cluster  size 
along  the  center  line  of  the  condensation  plume,  discussed  in  Sec.  4. 


Figure  12:  Comparison  of  Rayleigh  scattering  intensity  between  experimental  data  and  simulation 
results  along  the  plume  center  line.  The  dashed  line  indicates  intensity  of  a  flow  solution  without 
condensation. 


22 


Radial  Positic 


